//----------------------------------*-C++-*----------------------------------//
/*!
 * \file   State.hh
 * \author Jeremy Roberts
 * \date   05/23/2011
 * \brief  State class definition.
 * \note   Copyright (C) 2011 Jeremy Roberts. Based on Denovo.
 */
//---------------------------------------------------------------------------//
// $Rev::                                               $:Rev of last commit
// $Author::                                            $:Author of last commit
// $Date::                                              $:Date of last commit
//---------------------------------------------------------------------------//

#ifndef STATE_HH
#define STATE_HH

#include <vector>

#include "angle/AngleTraits.hh"
#include "angle/AngularMesh.hh"

#include "tran/SlabtranElement.hh"
#include "tran/BaseInput.hh"
#include "tran/CartMesh.hh"
#include "tran/Traits.hh"
#include "tran/Moments_Field.hh"
#include "tran/Psi_Cell_Field.hh"

#include "../utilities/DBC.hh"
#include "../utilities/SP.hh"
#include "../utilities/Vector_Lite.hh"

namespace slabtran
{

//===========================================================================//
/*!
 * \class State
 * \brief Problem state.
 *
 * The state of the problem at any point in the solution process can be
 * defined by the unknown flux moments (scalar flux and higher moments),
 * as these quantities are sufficient to describe reaction rates, which is
 * typically what we need (e.g. doses or fission rates).  For eigenvalue
 * problems, keff is also important.
 *
 * We also include the option to store the angular flux for each cell for
 * use in debugging and other applications.  This can be either the cell
 * centered angular or the cell faces angular flux, the latter of which
 * is required for CMFD.
 *
 * The moments are stored in \ref Moments_Field objects, with one such
 * object for each group stored in a vector.  The Moments_Field class
 * is templated on equation (for its number of cell unknowns) and
 * dimension (for dimension-specific indexing).
 *
 * The angular flux, when stored, is arranged similarly in a \ref
 * Psi_Cell_Field or Psi_Face_Field.  These are \em not yet
 * implemented.
 *
 * When needed, the underlying data array for a group can be accessed as
 * demonstrated by the following (I <b>think!</b>).
 * \code
 std::vector<Moments_Field<D> > group_moments;
 // build moments...
 double *localgrouparray;
 localgrouparray = &moments[g][0];
 * \endcode
 * This would be useful e.g. for filling a PETSc Vec object to use in one
 * of their Krylov schemes (rather than copying into another array).
 *
 * \todo Test whether referencing a Moments_Field object at [0] actually
 *       yields the array underneath.
 */
//===========================================================================//

template <class Equations, class D>
class State
{

public:

  //\{
  /// Useful typedefs for classes that use State
  typedef typename util::SP<State<Equations, D> >       SP_state;
  typedef BaseInput::SP_input                           SP_input;
  typedef typename Base_Mesh::SP_mesh                   SP_mesh;
  typedef typename AngularMesh<D>::SP_angularmesh       SP_angularmesh;
  // flux moments
  typedef Moments_Field<Equations, D>                   Moments_Field_t;
  typedef typename Moments_Field_t::SP_moments_field    SP_moments_field;
  typedef typename std::vector<SP_moments_field>        Scalar_Flux;
  // angular flux
  typedef Psi_Cell_Field<Equations, D>                  Psi_Cell_Field_t;
  typedef typename Psi_Cell_Field_t::SP_psi_cell_field  SP_psi_cell_field;
  typedef typename std::vector<SP_psi_cell_field>       Angular_Flux;
  //\}

  /*
   *  \brief Constructor.
   *
   *  \param    input           Input data.
   *  \param    mesh            Prebuilt mesh.
   *  \param    angularmesh     Prebuilt angular mesh.
   */
  State(SP_input input, SP_mesh mesh, SP_angularmesh angularmesh);

  /// \name Accessors
  //\{

  /*
   *  \brief Const accessor to a group moments field.
   *
   *  \param    g   Group of field requested.
   *  \return       Constant reference to a moments field.
   */
  const Moments_Field_t& get_phi(int g) const
  {
    Require(g >= 0 && g < d_number_groups);
    return *(d_scalar_flux[g]);
  }

  /*
   *  \brief Mutable accessor to a group moments field.
   *
   *  \param    g   Group of field requested.
   *  \return       Mutable reference to a moments field.
   */
  Moments_Field_t& get_phi(int g)
  {
    Require(g >= 0 && g < d_number_groups);
    return *(d_scalar_flux[g]);
  }

  /*
   *  \brief Const accessor to a group angular flux field.
   *
   *  \param    g   Group of field requested.
   *  \return       Constant reference to a psi cell field.
   */
  const Psi_Cell_Field_t& get_psi(int g) const
  {
    Require(g >= 0 && g < d_number_groups);
    Require(d_angular_flux[g]);
    return *(d_angular_flux[g]);
  }

  /*
   *  \brief Mutable accessor to a group angular flux field.
   *
   *  \param    g   Group of field requested.
   *  \return       Mutable reference to a psi cell field.
   */
  Psi_Cell_Field_t& get_psi(int g)
  {
    Require(g >= 0 && g < d_number_groups);
    Require(d_angular_flux[g]);
    return *(d_angular_flux[g]);
  }

  /*!
   *  \brief Return number of groups.
   */
  int number_groups() const
  {
    return d_number_groups;
  }

  /*!
   *  \brief Return size of flux moments for one group.
   */
  int moments_size() const
  {
    return d_moments_size;
  }

  /*!
   *  \brief Return size of angular flux for one group.
   */
  int psi_size() const
  {
    return d_psi_size;
  }

private:

  /// Input.
  SP_input d_input;

  /// Mesh.
  SP_mesh d_mesh;

  /// Angular mesh.
  SP_angularmesh d_angularmesh;

  /// Number of groups.
  const int d_number_groups;

  /// \name Fundamental unknowns
  //\{

  /// Multigroup scalar flux moments.
  Scalar_Flux d_scalar_flux;

  /// Multigroup angular flux.
  Angular_Flux d_angular_flux;

  /// Eigenvalue.
  double d_keff;

  //\}

  /// Size of one moments field
  int d_moments_size;

  /// Size of one angular flux field
  int d_psi_size;

};

} // end namespace slabtran


#endif // STATE_HH
//---------------------------------------------------------------------------//
//              end of State.hh
//---------------------------------------------------------------------------//
